Nonresonant Raman and inelastic X-ray scattering in the charge-density-wave phase 

of the spinless Fahcov-Kimball model 
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The dynamical mean-field theory formalism to describe nonresonant inelastic light and X-ray 
scattering in a charge-density-wave phase is developed and applied to the spinless Falicov-Kimball 
model on an infinite-dimensional hypercubic lattice at half filling. At zero temperature, the charge 
gap in the density of states is exactly equal to U ; increasing the temperature rapidly fills the gap 
with subgap states. The nonresonant response function for Raman and inelastic X-ray scattering 
shows peaks connected with transitions over the gap and transitions that involve subgap states; 
in addition, the spectra have significant changes in shape as the temperature is raised from zero 
to Tc. In the case of X-ray scattering (when both energy and momentum are transferred), the 
response function illustrates features of dynamical screening (vertex corrections) in the different 
(nonresonant) symmetry channels (j4ig and -Big); dynamical screening is also present in the Aig 
Raman signal. Finally, we derive and verify the first moment sum rules for the (nonresonant) Raman 
and inelastic X-ray response functions in the charge-density-wave phase and we discuss experimental 
implications for how the sum rules might be employed in data analysis. 

PACS numbers; 71.10.Fd, 71.45.Lr, 78.30.-j 



I. INTRODUCTION 

Charge-density-wave (CDW) systems possess a static 
rearrangement of the charge that is modulated by their 
ordering vector. Since the underlying ionic cores are 
charged, they will respond to this charge modulation 
from the electrons, and often create a distorted lattice 
structure that follows the modulated charge order of the 
electrons. This is often one of the easiest to measure 
signals of CDW order, namely the distortion of the unit 
cell due to the ionic displacement that goes hand-in-hand 
with the electronic charge modulation; it is more difficult 
to directly measure the electronic charge modulation in 
the material. 

In this work, we focus on signatures of the CDW or- 
der that are present in inelastic light scattering experi- 
ments on CDW systems. Since inelastic Raman scatter- 
ing is sensitive to different symmetry charge modulations 
(when polarizers are used on the incident and scattered 
light), it can provide information about the symmetry of 
the CDW state which is complementary to the results 
that would come from an elastic light scattering mea- 
surement such as optical reflectivity (which can measure 
only one symmetry). Similarly, because inelastic X-ray 
scattering also allows for an exchange of momentum by 
the scattered photon, we might anticipate interesting be- 
havior to occur when the ordering wave vector and the 
transferred momentum are the same. 

We develop all of the formalism to generalize the dy- 
namical mean-field theory (DMFT) approach to inelas- 
tic Raman and X-ray scattering in the situation when 
there is a CDW phase on a bipartite lattice with an or- 
dering wave vector equal to (tt, tt, . . . , tt); our formulas 
include all effects of vertex corrections including dynam- 
ical screening. While the formal development, in terms of 
the Green's functions, self-energies, and irreducible ver- 



tex functions, is completely general, and can be applied 
to any many-body model that has CDW order, such as 
the attractive Hubbard model or the Holstein model, we 
analyze the formalism for the specific case of the Falicov- 
Kimball model because the irreducible charge vertex is 
known exactly, and so we can provide an exact solution 
to the light scattering problem. In addition to deriving 
formulas for the light scattering spectra, we also exam- 
ine the first moment sum rules for these spectra, which 
are equal to expectation values related to the kinetic and 
potential energies of the material. These sum rules can 
be employed to aid in the data analysis of experiments, 
when higher-energy bands are well separated from the 
low energy band that undergoes the CDW order, as al- 
ready observed in systems that do not have CDW order, 
like in the normal state of SmBg at low temperature. 
They also provide an alternative way to directly measure 
the electronic order parameter of the CDW. 

We anticipate our results should be relevant to dif- 
ferent experimental systems that display charge-density- 
wave order via nesting on a bipartite lattice at half fill- 
ing, especially in compounds which are three-dimensional 
such as^ BaBiOa and Bai-ajK^jBiOa, because DMFT is 
most accurate in higher dimensional systems; it may also 
be relevant to some layered two-dimensional systems, at 
least in a semi-quantitative fashion. Our work also ex- 
tends recent results on transport and optical conductivity 
in CDW systems^"'^ to the realm of inelastic light scat- 
tering. Since inelastic light scattering experimental work 
on CDW systems has focused on Raman scattering of 
the soft phonon modes, the next step experimentally will 
likely be to examine the electronic scattering directly (ei- 
ther with Raman or with X-rays). Hence this work has 
the potential to be directly relevant to the next genera- 
tion of experiments in this area. 

The paper is organized as follows: in Sec. II, we derive 
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the formahsm for inelastic light scattering in a symmetry 
broken phase including explicit expressions for Raman 
scattering, inelastic X-ray scattering, and their first mo- 
ment sum rules; this formal development is appropriate 
for any many-body model of light scattering with local 
interactions. In Sec. Ill, we present our numerical re- 
sults for the example case of the spinless Falicov-Kimball 
model and discuss what signatures are likely to be seen 
in experiment. Our conclusions are presented in Sec. IV. 



local and nonlocal parts 

n 



(1) 



ijab 



where i and a = A or B are the site and sublat- 
tice indexes, respectively, and tfj is the hopping ma- 
trix, which is nonzero only between different sublattices 
(t^^ = tfj^ = 0). The local Hamiltonian is equal to 



II. FORMALISM 

Since CDW ordering is a static order, it is often well 

described by static models such as the Falicov-Kimball 
model.'* This model was introduced in 1969 to describe 
metal-insulator transitions in rare-earth compounds and 
transition-metal oxides. Since then, it has been studied 
widely within the DMFT community, primarily because 
it is one of the simplest many-body problems that ad- 
mits an exact solution^ (for a review see Ref. 6). The 
Falicov-Kimball model has two kinds of particles: mobile 
electrons and localized electrons. Mobile electrons hop 
from site to site with a hopping integral between nearest 
neighbors and they interact with the localized electrons 
when both sit on the same site (the interaction energy 
is U); we denote the mobile electron creation (annihi- 
lation) operator at site z by d| (rfj and the local elec- 
tron creation (annihilation) operator at site i by {li)- 
The model has commensurate CDW order at half filling 
and this is the main property we exploit here. Brandt 
and Mielsch worked out the formalism for calculating 
the ordered-phase Green's functions^ shortly after Met- 
zner and VoUhardt introduced the idea of the many-body 
problem simplifying in large dimensions.^ The CDW or- 
der parameter was shown to display anomalous behavior 
at weak coupling,^' and higher-period ordered phases 
have been examined on the Bethe lattice. Transport 
calculations in the commensurate CDW phase have also 
appeared recently.^'"' For concretcness, we will focus om* 
attention in the formalism section on the Falicov-Kimball 
model, but the light scattering formulas have a wider 
range of applicability. 



A. DMFT for the CDW ordered phase 

The hypercubic lattice is a bipartite lattice, implying 
that it separates into two sublattices (called A and B) 
with the hopping being nonzero only between the two 
sublattices. In this case, the model will display commen- 
surate (chessboard) CDW order when both the light and 
heavy particles are half-filled. This CDW order corre- 
sponds to the situation where the average filling of the 
electrons remains uniform on each sublattice, but changes 
from one sublattice to another (it is commensurate be- 
cause the lattice is bipartite here). We begin by writing 
the Falicov-Kimball model Hamiltonian as the sum of its 



(2) 



with the number operators of the mobile and localized 
electrons given by fiid = d^d^ and fiif = flfi, respec- 
tively. Note that we have introduced different chemical 
potentials for different sublattices. This is convenient for 
computations, because it allows us to work with a fixed 
order parameter, rather than iterating the DMFT equa- 
tions to determine the order parameter (which is subject 
to critical slowing down near Tc). Of course, the equi- 
librium solution occurs when the chemical potential is 
uniform throughout the system [ji"^ = ji^ and ji^ — iij ). 

We start with the definition of the lattice Green's func- 
tion 



Gt^ir) = -Tr[%e->'%,{T)d\{0)] jZ 
Z = Trexp[-/JW]. 



(3) 



Within a Feynman-diagram formalism, the Green's func- 
tion satisfies Dyson's equation (which in fact is a compact 
form of the diagrammatic series) 

Y}{^ + - ^tn^) + ttnCf^iuj) = 6ij6ab, (4) 



where w is a real frequency. The unperturbed band 
structure for the hypercubic lattice with nearest neighbor 
(NN) hopping satisfies 

D 

efe = - E ^t)^ exp[ik • (R,^ - Rjb)] = -2i ^ 

i — j OL=l 

(5) 

where Rj^ is a lattice vector for site i on sublattice A 
and a is the lattice spacing (we set a = 1). 

The first step of DMFT is to scale^ the hopping matrix 
element as t = t* j^^pD (we use = 1 as the unit of 
energy) and then take the limit of the infinite dimensions 
D — > oo. The self-energy is then local: 



S«j'(a;) = S«(a;)%^„6, 



(6) 



and in the case of two sublattices has two values 

and YP [uj). Now, we can write the solution of the Dyson 
equation (in a momentum representation) in a matrix 
form 



(7) 
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where z{uj) and the hopping term tk are represented by 
the following 2x2 matrices: 



tk 











(8) 



efc 
fife 



After substituting Eq. (8) into Eq. (7), wc obtain three 
equations for the different Green's function components 



(9) 
(10) 
(11) 



with Z{u!) defined by 

Zico) = ^[u; + S^(a;)][c. + - S^(a;)]. (12) 

These expressions agree with those of Brandt and 
Mielsch'^ even though our notation is somewhat differ- 
ent from theirs. 

The local Green's functions are now found to be 

G'^'^H = 4y Gr(c.) = ^^^±4f^^oo[^(-)], 



where 



dep{e) 



Z{uj) - e 



(13) 



(14) 



is the Hilbert transform of the nonintcracting density of 
states (DOS), which satisfies p(e) = exp(-e^/r^)/rV7r 
for the infinite-dimensional hypercubic lattice. 

The second step of DMFT is to map the lattice Green's 
function onto a local problem by means of the dynamical 
mean field. Since there are two sublattices, a dynamical 
mean field A"((jj) is introduced on each of them. As a re- 
sult, the local lattice Green's function on each sublattice 
becomes: 



G"°(a;) 



1 



(15) 



The third equation that closes the system of equations 
for G""(ti'), S"(a;) and X°-{uj) is obtained from the condi- 
tion that the local Green's function can be defined as the 
Green's function of an impurity with the same dynamical 
mean field X°'{uj). Such a problem can be exactly solved 
and the result is equal to 



G""(w) 



w + Md-A«(a;) w -|- /x^ - C/ - A»(a;) 



. (16) 



This last equation must be modified if one solves a differ- 
ent many-body model such as the Hubbard or Holstein 



model, as one needs to solve the relevant impurity prob- 
lem for the model being considered; the remainder of the 
algorithm is identical for other models. 

These equations are self-consistently solved numeri- 
cally. The iterative DMFT algorithm to calculate the 
lattice Green's function is as follows: we analytically 
continue all of the above formulas to the Matsubara 
frequency axis, because calculations along this axis are 
much more stable and converge faster than those on the 
real axis. Then, for a fixed value of the order param- 
eter An/ = — nj, one chooses and in such 
a way that + nj = 2n/ {nj- = 1/2 for half filling). 
With those fixed quantities, we now propose a guess for 
the self-energy on each sublattice, and then compute the 
local Green's fimction from Eqs. (12) and (13). Then 
we extract the dynamical mean field on each sublattice 
from Eq. (15), and find the local Green's function for the 
impurity from Eq. (16). This value is substituted into 
Eq. (15) to calculate the new self-energy. This proce- 
dure is repeated until the Green's function converges and 
we can calculate the filling of the conduction electrons; 
the chemical potential for the conduction electrons is ad- 
justed so that the average conduction electron filling is 
equal to one half. In order to find the correct equilibrium 
order parameter An/ at the given temperature, one next 
calculates the chemical potentials for the /-electrons on 
each sublattice via 



f,} = ---Tin— ^ 



(17) 
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T^ln 



1 - 



U 



iujn + IJ'd - A"(«a;„) 



where we introduce the fermionic Matsubara frequencies 

iujn = iTTT{2n+ 1). If these two chemical potentials are 
not equal, then the order parameter chosen initially is in- 
correct, and one needs to repeat the iterative loop with a 
new An/ to eventually satisfy the equilibrium condition 
where /U/ — /Z/ = 0; when this condition is satisfied, then 
An/ is the order parameter at that temperature. This al- 
gorithm, where the order parameter is fixed and we check 
for equilibrium by examining the chemical potentials, and 
then update the fixed order parameter to achieve the 
equilibrium solution, does not suffer from critical slow- 
ing down, which does occur if we instead fix the chemical 
potentials and iterate the equations until they converge. 
Generically, the DMFT equations can be solved with an 
order of magnitude less computer time than if we use this 
alternative approach. Finally, we repeat this iterative so- 
lution on the real axis, with the chemical potentials and 
fillings fixed at their now known values, which also is 
much more efficient than trying to do the entire calcu- 
lation on the real axis. For more complicated models, 
one most likely will need to fix the chemical potential 
and iterate the equations (which will be subject to criti- 
cal slowing down near T^), because one cannot solve the 
impurity problem with a fixed order parameter anymore. 
This does not create any serious problems, it just requires 
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more computer time. 

In Ref. 3, we already analyzed the evolution of the DOS 
in the CDW-ordered phase. We reiterate the main points 
which will be needed here. At T = 0, a real gap develops 
of magnitude U with square root singularities at the band 
edges. As the temperature increases, the system develops 
substantial subgap DOS which are thermally activated 
within the ordered phase. Plots of the DOS can be found 
in Ref. 3. Note that the singular behavior occurs for one 
of the "inner" band edges on each sublattice, and that the 
subgap states develop very rapidly as the temperature 
rises. 



B. Nonresonant inelastic scattering 

Now we develop the formalism for nonresonant light 
scattering in the CDW phase. We start from the stan- 
dard formula for the inelastic light scattering cross sec- 
tion derived by Shastry and Shraiman^^ 



contributions 



= {fhaAQ)\i) (22) 
\fMkf)\l) {l\3a{-ki)\i) 



+ 



{f\jai-ki)\l) {l\j0{kf)\i) 

ei-Si+ Wf 



with the sum / over intermediate states, and after substi- 
tuting into the cross section formula, one obtains three 
terms: a pure resonant term; a nonresonant term; and a 
mixed term (because it is constructed from the square of 
the scattering operator). 

The nonresonant contribution is 



(23) 



Eexp(-/3£i) . . 
5 liJlf^A^f - Si 



(18) 



The tilde denotes contractions with the polarization vec- 
tors: 



X \g{h)g{kf)el^eUf M'^0{q) 



1 = Y^(^alccA(l)(^0 

a/3 



(24) 



It describes the scattering of band electrons by photons 
with O = iVi — u)f and q = ki — kf being the trans- 
ferred energy and momentum, respectively, e*^^) is the 
polarization of the initial (final) states of the photons 
and £j(y^) denotes the electronic eigenstates. The quantity 
g{q) = {hc^ /VujqY^^ is called the "scattering strength" 
with Uq = c\q\. The scattering operator M{q) is con- 
structed from both the number current operator and the 
stress tensor which are equal to 



with the notation Oij = {i\0\ f) for the matrix elements 
of an operator O. (Resonant and mixed diagrams will be 
examined elsewhere.) 

The next step is to evaluate the summations in Eq. (23) 
via Green's function techniques. In general, such a proce- 
dure is nontrivial. But for the nonresonant contribution 
it is relatively straightforward^^' We start from the 
Matsubara function built on two time dependent stress- 
tensor operators 



M = E ^^iC^ + - 9/2) (19) X,,,(r - rO = TV [T.e-'^«7(r)7(r')] f^- (25) 



abk 



and 



7a,/3(9) = E 



abk 



dlik + g/2)4(fc - q/2), (20) 



respectively (in models with spin, and additional sum 
over the 2:-component of spin is required). Here tab{k) are 
the components of the 2x2 hopping matrix in Eq. (8). 
The interaction of an electronic system with a weak ex- 
ternal transverse electromagnetic field A is described by 
the Hamiltonian 



The imaginary time dependence of the stress-tensor op- 
erator is evolved (in the Hciscnberg representation) with 
respect to the equilibrium Hamiltonian because this is 
a linear-response calculation. The symbol %■ is a, time 
ordering operator. Further, we perform a Fourier trans- 
formation to the imaginary Matsubara frequencies. In 
thermal equilibrium, the two-particle correlation function 
depends only on the difference of the two time variables 
and our Matsubara frequency dependent function can be 
evaluated as 



^-t = -^E^'('^)-^(-^) 



(21) 



X7,7(«^) = E 
iJ 



exp(-/3£») 7»,/7/,i 
Z Sf — Ei — if 

X [1 - exp(/3(ei - £/))] . 



(26) 



J2M-khaAk + k')A0{-k'). 



kk' 



Performing an analytic continuation to the real axis iv 
± iO"*" produces the known expression 



The scattering operator M is then constructed from these 
interaction terms; it has both nonresonant and resonant 



2^g^{k,)gHkf) 



(27) 
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where we introduced the nonresonant response function pressed directly in terms of the generahzed polarizations 



TT 



(28) 



Now we have reduced the problem to that of finding 
the response function built on two stress-tensor opera- 
tors. Actually, such a function corresponds to a two- 
particle Green's function that will be shortly presented 
in Fcynman diagrammatic notation. The Fourier trans- 
form of the two stress-tensor correlation function can be 
represented as a sum over Matsubara frequencies of the 
"generalized polarizations" 



(29) 



where we introduced the shorthand notation lim.m+i = 
Hiiujmjii^m + ii^i) for the dependence on the fermionic 
iujm = iTrT{2m + 1) and bosonic ivi — i2-KTl Matsubara 
frequencies. In the case of the CDW ordered phase, the 
Feynman diagrams for the "generalized polarizations" 
'^m,m+i are shown in Fig. 1, where we introduce addi- 
tional sublattice indices a to L Here, we used the fact 



k-Cj/2, (ft), 



k-q/l. 1(0,. fiK-qn, »i>. 




k+q/2,ili3,„+iVi ^ I k'+q/2, Ka,,, +iv, 



FIG. 1: The Feynman diagrams for the generalized polariza- 
tions. Due to the properties of the irreducible charge vertex 
of the Falicov-Kimball model, we will have m — m' . 



that the total reducible charge vertex (shaded rectan- 
gle in Fig. 1) is a diagonal function of frequencies for 
the Falicov-Kimball model [see Eq. (42) below]; for other 
models, where the vertex is almost certainly no longer di- 
agonal, the analysis is somewhat more complicated. Now 
one can perform an analytic continuation to the real axis 
and replace the sum over Matsubara frequencies by an 
integral over the real axis: 



2m 



dujf{uj) 



(30) 



X n(a; — iO+, w + w;) — n(w -f iO"*", w -I- ii^; 
-I- Ii{uj - iiyi,uj~ iO+) - Il{uj + iiyi,uj + iO+) 



where /(w) = 1 /[exp(/3a;) -I- 1] is the Fermi distribution 
function. Then the nonresonant response function is ex- 



duj [f{uj) - fiiu + n)] (31) 



(27ri)2 

— OO 

Re|n(w - iO+, uj + n + iO+) 



The next step is to calculate these generalized polariza- 
tions. We consider both cases of inelastic light (Raman) 
and inelastic X-ray scattering. For Raman scattering, we 
can approximate q = because the optical photon wave- 
length is so large, whereas for inelastic X-ray scattering, 
the transferred momentum is nonzero q ^ 0- 



C. Raman scattering: q = 

The non-resonant Raman response function presented 
in terms of the generalized polarizations in Eq. (31) is 
reduced to the calculation of the Feynman diagrams in 
Fig. 1. As a result, our aim is to calculate the sum of the 
products of the one-particle Green functions calculated 
in DMFT and the charge vertices. Here, the momen- 
tum k enters not only through the band energy term 
Efc [see Eqs. (9-11)] but also through the stress-tensor 
factors, namely the derivatives d^e{k)/dkadkj3. Further- 
more, the stress-tensor operator is contracted with po- 
larization vectors e*'-'', see Eq. (24), which vary for the 
different symmetries. 

There are three symmetries often examined in experi- 
mental systems with cubic symmetry. The Aig symmetry 
has the full symmetry of the lattice and for the hyper- 
cubic lattice the incident and scattered light are both 
polarized along the same diagonal direction, so in large 
dimensions we take the initial and final polarizations to 
be = = (1, 1, 1,1,.. .). The stress-tensor amplitude 
in the case of Aig symmetry (for NN hopping) is equal 
to minus the band energy 



a/3 



^ a— 1 

(32) 



The Big symmetry is a d-wave-like symmetry that in- 
volves crossed polarizers along the diagonals. In this case, 
we take e* = (1, 1, 1,1,.. .) and = (—1, 1,-1, 1, . . .), 
so the stress-tensor amplitude is as follows 



fi\ ST' i f d'^^i.k) 



D 



a/3 



— ^(-l)"cosfc„ 

^ ^ a=l 

(33) 



Finally, the symmetry is another c?-wave-like symme- 
try rotated by 45 degrees; it requires the polarization vec- 
tors to satisfy e* = (1, 0, 1,0,.. .) and = (0, 1, 0, 1, . . .), 
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and for NN hopping there are no contributions to the 
nonresonant response in this channel. 

We start with the analysis of the Big symmetry, which 
is simplest case to examine. Here, the response function 
is determined only by the first term (bare loop) of the 
Feynman diagrams in Fig. 1 and there are no contribu- 
tions from the second one^^'^^ because the stress-tensor 
factor has momentum dependence that integrates to zero 
when multiplied by the local charge vertex and summed 
over all momentum. The expanded form of the diagrams 



and 



2^lk'^k,m'-^k,m+l — ^ Ik'^ k,m'-^ k,m+l ("^ ' ) 

k k 

1 Zjri+lFoo{Zm+l) ~ ZmFoc{Z,n) 

2 



72 _ 72 



Hence, the total expression for the generalized polariza- 
tion Ilm,m+i is 



k-q/2. ;co,. 



Foo{Z,n+l) FaoiZ„ 



abed jy 




H. 



72 _ 72 



(38) 



+ {iuj.m + Aid - S™)(iw„,. + ii^i + //f - ^^i+i) 

Zm+lFooiZyn+l) ~ ZyyiFoc{Z„i) I 



+ 2 



72 _ 72 



A B 

Y + Y| 

A B A ^ A 

k+q/2, iQ),.+iv, k+q/2, ia),„+(V, 

FIG. 2: Individual terms for the bare polarization in the or- 
dered phase. 

for the generalized polarization in the Big channel for 
the CDW chessboard phase is presented in Fig. 2 and is 
equal to 



Then, after substituting this expression for Tlni.m+i into 
Eq. (29) and replacing the summation over fermionic 
Matsubara frequencies by integrals over the real fre- 
quency axis, the total expression for the nonresonant re- 
sponse function equals 



— 00 

F4[Z(c. + 17)] F^[Z{cj)] 



(39) 



H™,™+; - — V 7fc ( G^-H Gffg 



N ^ 

k 



k—2-m k+^^m+L k—^^m k-\--^.7n 



X Re-! 



z*{Lu + n) 



+1 



(34) 



G 



G 



BA 



BA 



riBB nAA 



After substituting in the expressions for the Green's func- 
tion in Eqs. (9), (10), and (11), and the expressions for 
the 7fc amplitude from Eq. (33), the individual contribu- 
tions to Hm^,„_(_i at q = become 



X + Aif - S^(w)][c^ + r! + Aid - + ^)\ 

+ [cj + Aid - S-^(w)]k + r! + Aif - S^*(cj + ^)\ 

Z*[uj + n)F^ [Z{u + »)] - Z{w)F^ [Z{oj)] 
Z-^{lo)-[Z*{uj + VI)Y 
F^[Z{Lo + n)] F^[Z{lo)] 



"j\r 2-^'k'^k,m'^k,m+l — 2^'^'" ' /^d ^m) 



N Z-^ 

k 



(35) 



Z(^) 



Foci{Z, 



m+l j 



FoQ {Zyu) 



X {iuJm+m + fJ-d - ^m+l)- 



-im+l 



X! ^kGk!LGt^+i = 2 (*^'» + Md - ^ni) (36) 

k 

Foo{Zm+l) Foo{Zm) 



X (iuJra+m+l-id - ^m+l) 



Zm+l 



Z-n 
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z^{uj) - z^{uj + n) 

X (^[uj + Aif - S^(w)][tJ + + Aid - + ^)] 

+ [lu + i4 - s^(^^)]k + + Aif - S^(cj + n)] 

Z(uj + n)F^ [Zju; + O)] - Z{lu)F^ [Zjco)] \ 

z^{lj)- z^{uj + n) j 



One can check that this expression for the Raman re- 
sponse function (for Big symmetry) in the CDW phase 
is connected with the one for the optical conductivity^ 
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a a 



FIG. 3: The irreducible charge vertex becomes local in DMFT 
and, accordingly, the reducible charge vertex depends only on 
two sublattice indexes. 



by the Shastry-Shraiman relation^^ 



(40) 



indicating that this relation continues to hold even in 
the ordered phases. This formula holds for all models, 
because it does not depend on the vertex function. 

The case of Aig symmetry has both terms of the Feyn- 
man diagram of Fig. 2 contributing to the expression 
for the nonresonant response function. According to the 
form of the stress-tensor factor, the summation over mo- 
mentum of the bare loop yields 



n 



(1) 



1 I Zjn^lFoo^Z^n 



Z ni Fqo ( ■^m ) 



(41) 



m+l) 



which is different from the one for the Big symmetry in 
Eq. (38). 

The second term of the Feynman diagram in Fig. 1 de- 
scribes the charge screening effects through the reducible 



charge vertex, which is defined from the irreducible one. 
In the DMFT approach, the irreducible charge vertex Fq 
is local and different for different sublattices (see Fig. 3); 
nevertheless, it has the same functional form as in the 
uniform phase, which is equal to^^ 



m) = 



(42) 



"pa 



1 



for the Falicov-Kimball model (an explicit formula for 
other models is not known). This expression also fol- 
lows from the partially integrated Ward identity, derived 
by Janis^^. Accordingly, the reducible charge vertex in 
the CDW chessboard phase depends on two sublattice 
indexes and is defined by the Bethe-Salpeter equation in 
Fig. 4 



pa& _ X pa I T^pa \ ^ ac pcb 



(43) 



k-q/2, iCO™ 

a a a fl ^ c fy 



a a a a c b 

k+q/2, ![0„+iv, 



FIG. 4: The Bethe-Salpeter equation for the reducible charge 
vertex in the CDW chessboard phase. 



where we introduce the bare susceptibility 



at 



b /^ba 



(44) 



Now, the generalized polarization can be presented in a compact matrix form as follows 



n 



(2) 



^ ^ [ 7fc Ik] 



X T 



r-iAA 

r<BA 

k-^ 



yAA 

q^m.va+l 
fBA 



yAB 

q.ra.7n+l 
pBB 

q.ra.in+l 



AB 



AA 



riAB 

riBB 

K 2 ' 



BB 



BA 



(45) 



N ^ 

k' 



AA 

k'—%,m 
BA 



BA 



G 



AB 



BB ' 



k'~-^^Ta k' + -^^ni+l 



G 



AA 



-iAB 



Ik' 
Ik' 



The next step is to put q — Q, expand the expression via partial fractions with respect to the band energy and 
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calculate the sums over momentum k. After some tedious algebra, we obtain the final expression for n^^^_|_; 



n 



(2) 



1 [Zm+lFoo{Zm+l) — ZmFao{Z„ 



m,m+l A 



[Zm - Z^+l] - («Wm + - + Wl + Hd - ^m+l) 



-'m+l 



- Gil, 



+ [i(2a;„ + v{) + - S;^ - S;^^,] ^^^^ - Z^F^{Z^)\ 

yB _ yB / 



Foo{Zm+l) Fco{Zm) 



-'m+l 



yA _ yA 
nAA _ r<AA 



+ [i(2a;„ + v{) + 2^1 - - E^+;] ^£ - Z^F^{Z^)\ j I 



(46) 



where 



D. X-ray scattering: q 



Am,m+J = det ||i5afc 



(47) 



comes from the solution of the Bethe-Salpctcr equation 
in Eq. (43). Finally, the total expression for the gen- 
eralized polarization is obtained as the sum of the two 
contributions 



n. 



= n; 



(1) 



n 



(2) 



(48) 



Next, we perform an analytical continuation (which is 
straightforward because we have the appropriate func- 
tional forms which allow us to replace Matsubara fre- 
quencies by real frequencies) and substitute into Eq. (31) 
which yields the final expression for the nonresonant Ra- 
man response function in the ^ig channel. This step 
is completely straightforward, so we do not write down 
the final expressions in terms of integrals over the real 
frequency. 

For other many-body models, such as the Hubbard, 
or Holstein models, the analysis is more complicated be- 
cause the vertex is not diagonal, and the analytic continu- 
ation will not be possible on a formal level, if the charge 
vertex cannot be expressed as functions of the Green's 
function and self-energy. Nevertheless, one can perform 
the analytic continuation numerically, and thereby solve 
the problem, or one can approximate the vertex using dif- 
ferent analytic approximations (such as those from per- 
turbation theory) and then formally complete the ana- 
lytic continuation. 



In the case of inelastic X-ray scattering, the incident 
photon exchanges both energy and momentum with the 
electronic matter. The entire formalism derived for Ra- 
man scattering remains the same as described above and 
there is no need to rewrite it for this case. The only dif- 
ference is in the summations over momentum. The Feyn- 
man diagrams in Fig. 1 together with the Bethe-Salpeter 
equation in Fig. 4 contain several momentum summa- 
tions which can be evaluated separately^®. First, the 
bare susceptibility in Eq. (44). which enters the Bethe- 
Salpeter equation for the total charge vertex in Eq. (43), 
contains the following components 



2ZmZm+l 



X [Xo{Zm, Zm+l, q) - XoiZm, -Zm+l, Q)] 



(49) 



and 



"^ZmZm+l 

[xo{Zm, Zm+l, q) — Xo{Zm, —Zm+l, Q)] 



y^S ^ BA 

= \ [Xo{Zm, Zm+l, q) + Xo{Zm, -Zm+l, ?)] , 



(50) 



(51) 
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where 



Xo{Zm, Zm+l, q) = Xo{-Zm, -Zm+l, Q) (52) 



Zm - <^u_l Zm+l - <^i,.a. 



+00 



y/T^^ J Zm+l - e 



k+7^ 

de , ( Zm — 



Here, the function Fac{Z) is the Hilbert transform of the 
hypercubic density of states as defined in Eq. (14) and all 
the transferred momentum dependence is only through 
the quantity 



1 ^ 



(53) 



p=i 



The second diagram in Fig. 1 contains summations 
over k and k' which involve stress-tensor amplitudes 7fe 



x; 



(54) 



aA 



GJia r^i. 

k—2,m k+2 ,Tn+l fe— k+2,m+l 



and there are two different terms: 



Xq,m,m-\-l n ry 



B v-B 



(55) 



X [x!o{Zm,Zm+hq) + X'o{ZmT-Zm+l,q) 
B 

m+l 



+ 



iuJm+l +IJ-d 



'^Zm+l 

X [x'oiZm, Zm+l, q) - x'oiZm, -Zm+l, ?)] 



and 



(56) 



^q,m,m+l r\ ry 

X \)<lQ{Zm,Zm+uq)^xlQ{Zm,-Zm+l,qf\ 

^ iuJm+l + \J^i-^t+l 
'^Zjn+l 

X [x'oiZm, Zm+l, q) - x!o{Zm, -Zm+l, Q)] , 



where 



N 
X' 



e^_| Zm+l - Efe+I 



= _^ ^ I [Zm + Zm+l] Xo{Zm,Zm+l,q) 
+ F^[Z{iWm)] + F^[Z{iuJm+l)]]- 

Here the new momentum dependent quantity X' is 



cos 



(58) 



p=i 



with Up = 1 for Aig symmetry and Up = (—1)'' for _Big 
symmetry. Now we can find exact expression for the 
vertex corrections defined by Eq. (45) in the following 
form 



H 



(2) _ _ 



1 



A rpy.A AB rpyB B 

A.Q,m,m-\-l m,m+/ Aq,m,m+/ m.m-\-lAq,m,m-\-l 



(59) 



\^ ^ m,m+lXq^m,m+l) ^ 



' .m^7n+l J ^ ^ m,m+lXq.m.m+l 
Aq,m,m+l V m,m+lA.q,m,m+l J m,m+lA.q,m,m+l 



,B 



-L TF^ v^^ TF^ v"^ 

A.q,m,m+Z m,m+l'^q,m,m+l m,m+lA.q,m^m+l 



where 

Aq,TO,m+; (60) 
= h - TF^ v^^ h - TV^ v-^^ 

V m,m+lA-q,m,m+lJ \ m,m+lA-q,m,m+l/ 

- TF^ v^^ TF-^ v^^ 

m,m+lA-q,mym+l m,m+l^q,m,m+l 

Finally, the bare loop contribution of the first diagram 
in Fig. 1 contains summations over momentum A; of the 
product of two Green functions and the square of the 
stress-tensor factor. It is equal to 



H 



-I- 



(1) 

q,m,ni-\-l 

1 



= XoiZm, Zm+l, q) + XoiZm, -Zm+l, ?) (61) 



"^ZmZm+l 

+ (iiOm +l^d - ^m){i^m+l + l^d - ^m+l) 
X [Xo{Zm, Zm+l, q) - XoiZm, -Zm+l, ?)] 

and expressed in terms of xo as follows 

XoiZm, Zm+l, q) = Xoi—Zm, —Zm+l, q) 



(62) 



-y 



ll- 



^ Ji*' t*^X'^ , [Zm + Zm+lYX'^ 

XoK^m, ^m+l,q> 



l + X 



+ 



(1+X)2 
X'2 



\_Zm + Zm+l ] {F^[Z{iuim)]+F^[Z{iu;m+l)]] 



x'oiZm, Zm+l, q) = -x!oi-Zm, -Zm+l, Q) (57) l + X 



+ — [ZiiUm)] + [ZiitOm+l)] - 2} . 



The expressions for x'o and Xo derived above appear to be 
different from the ones given in Ref. 18. In fact, they are 
identical (but require some significant algebra to show 
this); the forms presented above are more convenient for 
numerical calculations. 

In contrast to Big Raman scattering at q = which is 
determined only by the bare loop contributions (Fig. 2), 
in the case of inelastic X-ray scattering, we have both 
terms contributing for all symmetry channels. The dif- 
ferent symmetry channels are distinguished only by the 
different X' factors, and, as a result, different x'o and 
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Xq functions. All further numerical calculations arc per- 
formed by exploiting these three quantities, but the total 
scheme remains the same. As a result, the generalized 
polarization in Eq. (59) is described in terms of the xo, 
Xo and Xo functions and applying further analytic con- 
tinuation to the real axis one can obtain the nonresonant 
inelastic X-ray scattering response functions. The final 
expressions are too long to be presented here. 



E. Nonresonant inelcistic X-ray scattering sum rule 

The sum rule for the nonresonant inelastic scattering 
response function is as follows: ^^'^'^ 

+ CXD 

1= j dnnxN{n) = ^{[^\q)[H,^{q)]]), (63) 



where for the case of CDW ordering 

7(9) = E E ^^^e^«(«?-«?)e-^l W+«?)4'i.., (64) 

ab ij 
ab ij 

and the momentum Q determines the symmetry channels 



fies 



Iu= [ diof{iv)lmlY.[E''{u;)-Un'}] 



{1-X")G'"'{lo)+X'^X''{lo) 



+ X'^ [Zioj)F^[Zioj)] - 1] [E^(w) - S^(a;)]' 



(68) 



-Ut*-^{l-X''){nf-nf){ni-nn 



(This last term is model dependent and would be dif- 
ferent for the Hubbard or Holstein model. We do not 
provide those formulas here.) The first contribution in 
braces has the same shape as the potential energy contri- 
bution of the sum rule in the uniform phase. The other 
terms appear only in the CDW phase and are propor- 
tional to the square of the CDW order parameter (An/)^. 

By examining different points in the Brillouin zone 
(BZ), one can extract information regarding the potential 
and kinetic-energy contributions or of the order parame- 
ter. For instance, in the case of Raman scattering (q = 0, 
X = 1) we have contributions only from the potential- 
energy term {Ix = 0), which are different for the A^g 
symmetry {X' = 1) 



Q = 



(7r,0,7r,0,...) 



for Aig 
for Big 



(65) 



After calculating all required commutators, taking the 
large dimensional limit, and performing some cumber- 
some transformations (see the Appendix), we obtain a 
sum rule (first moment of the nonresonant inelastic X-ray 
scattering response function) which contains two contri- 
butions 

I = Ik+Iti. (66) 
The first contribution comes from the kinetic energy term 

+ CXD ^ 

Ik = '2{1-X)J diofiiv) Im I ^ [Z{uj)Fo, [Z{u;)] - l] 

— OO 

(67) 



-X' 



Z\u;) 



[Z{w)F^[Z{u;)]-l]+^ 



and is similar to the one in the uniform case.^^ The other 
one originates from the potential energy term and satis- 



-l-oo ^ 

/n= f (Ljf{u))lml^[j:''{Lo)-Un'}]X''{Lo) (69) 



+ [Z{u>)F^[Z{u)] - 1] [S^(a;) - S^(a;)]- 



and for the Big symmetry (X' = 0) 



+00 ✓ 
Iu= J rfa;/(a;)Im|^ [5]«(a;) - Unj] ^G«"(a;)| 



(70) 



^Ue^nf-nf){ni-n^). 



For other points in the BZ (inelastic X-ray scattering), we 
have contributions from both the kinetic and potential- 
energy terms. For instance, for the case of Big symmetry 
along the BZ-diagonal [q = (q, q, q.q, . . .), — 1 < A' < 1, 
X' = 0] and for all symmetry channels at the BZ corner 



11 



' TT TT TT TT 

^ a ' a ' a ' a ' " 



.), X = -1, X' = 0] we have: 



+ CXD 



Ik ^ 2(1 -X) I dLofiu;)lm{ — [Z{Lo)Fo,[Z{Lu)]-l\ 



t*2 



(71) 



+ 00 



(72) 



One can see, that in this case the kinetic-energy contri- 
bution is equal (up to an overaU constant) to the average 
kinetic energy which also enters the sum rule for optical 
conductivity.'^ 



III. NUMERICAL RESULTS 

We begin with an analysis of nonresonant Raman scat- 
tering in the CDW phase. We present results for the 
cases of a dirty metal with U — 0.5 (Fig. 5), a near criti- 
cal Mott insulator with U = 1.5 (Fig. 6), and a moderate 
gap Mott insulator with U = 2.5 (Fig. 7). 

In Fig. 5, we plot the Raman response function for 
different temperatures in the case of a dirty metal with 
U — 0.5. At temperatures higher then the critical one 
for CDW order, we see the expected behavior for a dirty 
metal: namely, there is a peak at low energy and a spread 
on the order of the metallic bandwidth. The system does 
not have a low energy Fermi-liquid peak, because it is 
not a Fermi liquid. Below the critical temperature, when 
the CDW gap arises, the shape of the response function 
changes significantly. The main peak is shifted to higher 
frequency at ~ J7, which corresponds to transitions 
between the lowest band at lu < —U/2 and the upper 
band at > U/2 (see the DOS in Ref. 3). Two addi- 
tional peaks at lower frequencies correspond to the tran- 
sitions from the upper and lower bands onto the subgap 
states and between the subgap states (which are present 
for a wide range of temperatures below Tc but above 
T = 0). Because the subgap DOS vanishes at T = 0, 
these peaks must vanish with T ^ 0. In addition, be- 
cause the self-energy becomes a frequency independent 
constant on each sublattice at T = (0 on one sub- 
lattice and U on the other), the irreducible charge ver- 
tex, and hence the vertex corrections, also vanish at low 
temperature. In panel (a), we plot the nonresonant re- 
sponse function for the Big symmetry. In this symmetry 
channel, there is a sharp main peak with a square root 
singularity at T = 0. This behavior was already seen 
in the optical conduct ivity'^, and follows for the Raman 
scattering directly from the Shastry-Shraiman relation in 
Eq. (40). For the Aig channel, as plotted in panel (b), 
the response is much smaller and smooth (without sharp 
singularities) and there are two reasons for this. At high 
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FIG. 5: (Color online) Nonresonant Raman response function 
for the two symmetry channels [(a) being the Big channel and 
(b) being the Aig channel] in a dirty metal with U = 0.5. The 
set of curves corresponds to a range of temperatures from the 
uniform to the ordered phase {Tc ~ 0.034). 



temperatures, in the uniform phase and just below Tc, we 
have the effects of dynamical charge screening for the Aig 
scattering channel which suppresses the total response. 
On the other hand, the charge vertex is proportional to 
U'^n'j{l — vfj) (see Ref. 16) and in the charge-ordered 
phase, where ^ 1 and — > for T ^ 0, it decreases 

rapidly as e"^"/-^, where i?a is the activation energy re- 
quired to flip the occupation of the /-state at a single 
site. This thermally activated disordering of the chess- 
board phase also gives rise to the subgap states. As a re- 
sult, the vertex contributions (dynamical charge screen- 
ing) become negligible at low temperatures [for example, 
increasing the response in Fig. 5(b)] and the total Raman 
response is determined by the bare loop contributions 
only. The expressions for the bare loop contributions 
[Eq. (38) for Big symmetry and Eq. (41) for Aig sym- 
metry] follow from Eqs. (61) and (62) by choosing q ~ 
{X = \) and AT' = for the Big symmetry and X' — \ for 
the Aig symmetry; the terms with AT'^ in Eq. (62) deter- 
mine the difference of the bare responses in the different 
symmetry channels. From a mathematical standpoint, 
the presence of the additional terms due to a nonvanish- 
ing X' , removes the singularity at the CDW gap edge 
which is in xo- From a physical standpoint, the different 
symmetries respond differently to the charge excitations. 
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FIG. 6: (Color online) Nonresonant Raman response func- 
tion for the two symmetry channels in a near critical Mott 
insulator with U = 1.5. The set of curves corresponds to a 
range of temperatures from the uniform to the ordered phase 
(Tc ^ 0.075). 
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FIG. 7: (Color online) Nonresonant Raman response func- 
tion for the two symmetry channels in a moderate gap Mott 
insulator with U — 2.5. The set of curves corresponds to a 
range of temperatures from the uniform to the ordered phase 
(Tc ^ 0.072). 



The Aig channel responds with the full symmetry of the 
lattice, and the effect of the modulated CDW tends to 
average out this response so that the singular feature dis- 
appears, while this does not occur for the Big channel, 
where a greatly enhanced response occurs near the charge 
gap. Hence, even though the vertex corrections vanish, 
which normally are required to guarantee that the system 
can screen long wavelength uniform charge fluctuations, 
the averaging effect of the CDW guarantees that the sys- 
tem can continue to screen these charge fluctuations even 
when the vertex corrections vanish. For example, this is 
precisely how the uniform charge susceptibility will van- 
ish in the CDW phase at T = 0, which is required by 
the equations of motion for the total charge, and occurs 
due to the vertex corrections in the normal state above 
Tc, below Tc, since the vertex corrections are suppressed, 
this averaging feature takes over and allows the system 
to continue to effectively screen out those charge excita- 
tions. 

In Fig. 6, we plot the results for a near-critical Mott 
insulator with U ~ 1.5. The basic results remain quite 
similar to the metallic case. We see the response function 
change dramatically as the system orders, with complex 
behavior at low temperature and low energy due to the 
subgap states, and then finally leading to the square root 



singularity in the Big channel and smoother behavior in 
the Aig channel, with no singularity, and significantly 
reduced spectral weight. The main change is the energy 
scale, since the gap is always identically equal to C/ at T = 
0, and this is reflected in the "pushing" of the spectra 
to the right. As we go from a near-gap insulator to a 
moderate-gap insulator with U = 2.5 (Fig. 7), we once 
again see similar kinds of behavior. In particular, we 
observe three peaks: the main CDW-gap peak at Q = U 
is sharp for the Big symmetry in panel (a) and smoothed 
for the Aig symmetry in panel (6) and the two low-energy 
peaks have strong temperature dependence. 

For nonresonant inelastic X-ray scattering, we investi- 
gate the behavior of the response functions for the differ- 
ent transferred momentum values q in the first Brillouin 
zone (BZ). Because all the momentum dependence en- 
ters only through the parameters X and X', we must 
first understand their behavior in the BZ. We want our 
results to make contact with real physical systems, like 
a two-dimensional system, so we choose the following 
paths in the first BZ: the zone diagonal (zd) path lies 
in the so-called E-direction with q = {q, g, q,q, . . .) and 
— 1 < X < 1; the zone edge (ze) path lies in the Z- 
direction with q = {- ,q, -,q, . . .) and —1<X<0, 
and then continues along the zone edge path in the A- 
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direction with q = {q, 0, q, 0, . 
results are depicted in Fig. i 



and < X < 1. These 
The corresponding de- 
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Z: q=(TT/a,q,TT/a,q, 



FIG. 8: Schematic of the first Brillouin zone with the high 
symmetry points labeled. Although we work in infinite dimen- 
sions, we are trying to make contact with the two-dimensional 
BZ. 
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FIG. 9: Plot of X and X' along the zone diagonal path and 
zone edge path in the first Brillouin zone. 

pendence of X and X' along these paths is plotted in 
Fig. 9. One can see, that along the Z-dircction, the X' 
value and, as a result, the response functions, are the 
same in both symmetry channels. For other directions, 
they are different. In addition, X' — Q along the zone 
diagonal S-direction for the Big symmetry and the cor- 
responding response function is determined only by the 
bare contributions with no vertex corrections (or dynam- 
ical screening) entering. 

Having established the values of X and X' that we are 
using, we now show our numerical calculations of the non- 
resonant inelastic X-ray response functions for the case 
of a dirty metal with U — 0.5 at different temperatures 
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FIG. 10: Nonresonant X-ray scattering response function in 
the -Big channel for U = 0.5 along the zone diagonal and zone 
edge of the first Brillouin zone. The set of curves correspond 
to temperatures T = 0.04, 0.03, 0.02, and 0.015. 
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FIG. 11: Nonresonant X-ray scattering response function in 
the A\g channel for U = 0.5 along the zone diagonal and zone 
edge of the first Brillouin zone. The set of curves correspond 
to temperatures T = 0.04, 0.03, 0.02, and 0.015. 



and transferred momentum. In Fig. 10, we present re- 
sults for the -Big symmetry and in Fig. 11 for the Aig 
symmetry. At the zone center {X = 1), the response 
is the Raman scattering (see figures above) with sharp 
features in the -Big channel and with a strong suppres- 
sion in the Aig channel. When we move away of the 
zone center, first of all, the sharp square root singular- 
ity at = [/ in the .Big channel is rapidly replaced by 
a step-like response with a strong enhancement at the 
Brillouin zone corner AT = — 1, when the transferred mo- 
mentum coincides with the CDW wave vector and we 
have effects of nesting present. For the Aig symmetry, 
we have a different scenario: there is a continuous en- 
hancement without any sharp features, when we move 
along the zone diagonal and there is a continuous devel- 
opment of a step-like feature, when we move along the 
zone edge with a strong enhancement at the zone corner 
also due to nesting. In addition, the screening due to 
the vertex corrections vanishes there for all temperatures 
(since the -Big and Aig response functions are identical 
and have no vertex corrections there [X' = 0]). In both 
cases, there is a large enhancement of the scattering re- 
sponse function as we move from the zone center to the 
zone corner; this occurs because the system, as a whole, 
is much more effective at screening out uniform charge 
fluctuations than those modulated in space. This shows, 
in particular, that scattering at the ordering wave vector 
is enhanced in the CDW system. 

Because the results for the near critical Mott insulator 
with [/ = 1.5 are similar to the results for the other two 
U values, we do not show them here. But, we do plot the 
results for a small gap Mott insulator, with U — 2.5 in 
Figs. 12 and 13). Here, we continue to see similar behav- 
ior to what is seen for U = 0.5, namely, the character of 
the response changes rapidly as we move away from the 
zone-center, the differentiation of the results for different 
symmetry channels is reduced, and the results coincide 
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FIG. 12: Nonresonant inelastic X-ray scattering response 
function in the Big channel for U = 2.5 along the zone di- 
agonal and the zone edge of the first BZ. The set of curves 
corresponds to the temperatures T = 0.08, 0.07, 0.06, and 
0.04. 
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FIG. 14: Big sum rules as a function of temperature for 
U = 2.5. (a) The kinetic-energy contribution Ik/{1 — X): 
1 — zone-diagonal S-direction (—1 < X < 1, X' = 0); 2 — 
X-point at zone-edge {X — 0). (b) The potential-energy con- 
tribution In: 1 — zone-diagonal E-direction ( — 1 < X < 1, 
X' = 0); 2 — zone-edge X = —0.5; 3 — zone-edge X = 
(X-point); 4 — zone-edge X — 0.5. The thin lines correspond 
to the uniform solution artificially continued below Tc- 



FIG. 13: Nonresonant inelastic X-ray scattering response 
function in the Aig channel for U — 2.5 along the zone di- 
agonal and the zone edge of the first BZ. The set of curves 
corresponds to the temperatures T = 0.08, 0.07, 0.06, and 
0.04. 



at the zone corner. We also see an enhancement of the 
signal and a generic broadening of the peaks as we move 
from the zone center to the zone corner. 

Since we have derived first-moment sum rules for all of 
the response functions, we checked our numerical results 
by integrating the first moment of the response function 
and comparing that answer to the results of the moment 
sum rule expectation values, which are evaluated on the 
imaginary axis. In all cases we examined, we achieved 
essentially perfect agreement, with errors less than 0.1%, 
and arising primarily from the discretization we used in 
our frequency grid for the numerical integrations. 

But the sum rules can actually tell us more about the 
system. One of the hallmarks of the /-sum rule for the 
optical conductivity is that the sum rule is fixed and does 



not change with temperature or interaction strength, so 
spectral weight is never lost or gained. In a projected 
low-energy model, this result no longer holds, and the 
low-energy spectral weight can change with temperature 
or U, but, as is often the case, the changes are quite 
small at low temperature. We can of course investigate 
this for our system in the CDW phase, by examining how 
the sum rule evolves for different parameters. We begin 
with a plot of the sum rule for the case of strongly cor- 
related insulator U = 2.5 in the Big channel in Fig. 14 
and for the Aig channel in Fig. 15. One can see, that 
for such values of Coulomb interaction the main contri- 
bution to the sum rule comes from the potential-energy 
part. The momentum dependence of the sum rule in the 
Big channel is weak for the potential-energy contribution 
and strong for the kinetic-energy one [notice the 1 — X 
factor in Eq. (67)]. For the Aig channel, both contri- 
butions have strong momentum dependence. For both 
symmetries, the largest values of the sum rule (total and 
for each contribution) are observed at the BZ corner M- 
point {X = —1) in both the uniform phase and the CDW 
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FIG. 15: Aig sum rules as a function of temperature for U — 
2.5. (a) The kinetic-energy contribution Ik /(I — X): 1 — 
zone-diagonal X — 0; 2 — M-point at zone-corner {X = — 1); 

3 — X-point at zone-edge {X = 0). (b) The potential-energy 
contribution In'- 1 — F-point at BZ center {X = 1); 2 — 
zone-diagonal X — 0; 3 — BZ corner X = —1 (M-point); 

4 — zone-edge X — (X-point). The thin lines correspond 
to the uniform solution artificially continued below Tc. 



phase, as could have been guessed due to the enhance- 
ment of the overall spectral functions we observed above 
(once again, in the CDW phase, we see an additional 
enhancement due to the ordering). The increase in the 
sum rule below Tc is linear in Tc — T and proportional 
to the square of the CDW order parameter (An/)^; this 
imphes that if an experimental system has a nice sepa- 
ration between the low and high energy bands, then one 
could use this spectral weight to measure the order pa- 
rameter as a function of temperature. For small values of 
U (see Figs. 16 and 17), the kinetic-energy contribution 
gives the main contribution into the total sum rule. The 
kinetic-energy contribution continues to display strong 
momentum dependence and for some momentum its tem- 
perature dependence becomes quite nonlinear below Tc, 
as we already saw for the optical sum rule.^ 



FIG. 16: -Big sum rules as a function of temperature for U 
0.5 (we plot the same cases as in Fig. 14). 



IV. CONCLUSIONS 

In this work, we developed the formalism (within 
DMFT) to calculate the nonresonant inelastic Raman 
and X-ray scattering in the case when the system de- 
velops CDW order at low temperature. The formalism 
is a straightforward generalization of the results in the 
paramagnetic phase, but requires a careful accounting of 
the different sublattices and how they enter into the di- 
agrammatic expansions, and hence is technically quite 
challenging. We also derived first-moment sum rules for 
these spectra and related the sum rules to different ex- 
pectation values that can be immediately calculated. We 
find that the sum rules relate to the potential energy in 
some cases, while in other cases, both the kinetic energy 
and the potential energy terms enter into the expectation 
values (and also the CDW order parameter). 

We applied our formalism to the case of the spin- 
less Falicov-Kimball model because the charge vertex is 
known exactly for that system, and hence we can find an 
exact numerical solution to the light scattering response 
functions. The main numerical result that we find is that 
there is very strong temperature dependence that sets in 
once we pass through Tc- This occurs because the sys- 
tem rapidly depletes subgap states as it forms the CDW 
gap, and then develops a square-root singularity due to 
the pile-up of states at T = 0. These features can be im- 
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FIG. 17: Aig sum rules as a function of temperature for U ■ 
0.5 (we plot the same cases as in Fig. 15). 



mediately seen in the light scattering response functions, 
but are symmetry selective. When vertex corrections act 
to screen the light scattering at high temperatures (near 
and above T^), the square root singularity is suppressed 
in the Aig channel, as is the overall magnitude of the 
light scattering signal. At low temperatures, the most 
important effects are due to nesting of the transferred 
momentum in BZ and due to an effective screening of uni- 
form charge fluctuations which arises due to an averaging 
effect over the modulated charge distribution of the sys- 
tem. The qualitative shape of the response function for 
inelastic X-ray scattering, where the spectra has almost a 
discontinuous jump near the gap, is an unexpected result, 
that occurs when one combines the square-root singular- 
ity associated with the pile-up of the density of states 
near the gap edge with the nontrivial nesting effects and 
the dynamical charge screening effects of the many-body 
system. While we see an enhancement of the response, a 
broadening of the spectra, and an increase in the magni- 
tude of the sum rule as we move from the zone center to 
the zone corner, we do not see any dramatic changes in 
the shape of the spectra associated with the fact that we 
can transfer momenta that is equal to the ordering wave 
vector of the CDW but we do see a significant overall en- 
hancement of the signal. This turns out to be similar to 
what was seen in the dynamical charge susceptibility of 
the model as one approaches Tc from above^^, and may 
be related to the fact that the Falicov-Kimball model has 



a reducible charge vertex that assumes very different be- 
havior for dynamical charge fluctuations as it does for 
static charge fluctuations, which give rise to the underly- 
ing CDW order. If true, then we would anticipate even 
larger effects in models where the charge vertex is not 
decoupled in this fashion, such as the Hubbard or Hol- 
stein model, but resolving this question is a problem for 
the future. 

Our numerical work focused on the case of half filling. 
One might ask what would happen away from half-filling. 
While it is true that the CDW phase can be seen as the 
first ordered phase as we go from the normal state to an 
ordered state at Tc,^^ we do not know whether the AB 
ordered phase survives all the way to T = or whether 
there are subsequent phase transitions, perhaps to in- 
commensurate phases as T is further reduced. For this 
reason, we have not chosen to solve such problems in this 
work. In the high-temperature phases, where the system 
is ordered in the AB CDW, the chemical potential would 
need to be located outside of the gap, and so we would 
expect to see more response at low energies, but as the 
T was further lowered, we expect incommensurate order 
to enter, and for the system to have a well developed 
gap, so that the results would most likely look similar to 
those shown here. On the other hand, there is another 
possible scenario at low temperatures for some densities 
of the mobile and localized electrons when, instead of the 
incommensurate order, the phase separation into chess- 
board and uniform phases can take place. ^^'^"^ In this case 
the total response will be a sum of the responses for the 
chess-board and uniform phases weighted by the volume 
fractions of these phases. 

We believe our results will be most relevant to elec- 
tronic Raman or X-ray scattering on CDW ordered sys- 
tems in three dimensions. So far, most of the Ra- 
man scattering work has focused on understanding how 
phonons behave as one passes through the transition, 
including the behavior of the phonon softening for the 
CDW mode^"*. We hope that our results will inspire ex- 
perimental groups to also consider examining electronic 
Raman scattering in CDW systems to see whether they 
display the kinds of features that we showed here. 

In the future, we will generalize the resonant light scat- 
tering formalism to the CDW phase and examine what 
modifications enter into the response functions in that 
case. 
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APPENDIX: SUM RULE DERIVATION 



Here we introduce the notations A = B and B = A and 
use the fact that the hopping integral connects only sites 
which belong to different sublattices. 



In this appendix, we present details for the derivation 
of the first-moment sum rules of inelastic light and X-ray 
scattering in the ordered CDW phase. To begin, we must 
evaluate the first commutator in Eq. (63) which yields 

[H,l{q)] = i7^t«/e^'3(«?-«J)e-'f(«?+«|) (A.l) 

ija 

ijla 

_ iQ{Rt-Rf) -i^iRt+Rf)] it i. 



The second commutator now gives 



ijla 



^ ' [riij.- — n 



/ - "If) 



dldja (A.2) 



j.aa 



iQ{Ri-Rf+R^-Rl)^ii(Rt+Rf-Rl-Rl) _ piQiR'^-R^) pi^iRf-R^) 



Ej.aaj.aaj.a 

ijlna 

. g«Q(K--Rf+K--Rj)e-*§«+«f-K--R|) _ g^QiRf-Ri)(,-^i{Rf-Rl)] dl^dja- 



Next, we use the fact that the hopping is allowed only between NN sites, and we replace j = i + 5 in tfj, where S runs 
over all of the NNs of site i, to obtain 



„-i%(d+S') (f^a f^a \ i!l{S+S') (-a ~a \ 



(A.3) 



-iQ{S'+S") 



iSS'5"a 



di^adi+S+S'+S",a 



The first term contains expectation values of three operator products which can be calculated by introducing an 
auxiliary field /Uj ^ /Xj -|- j at site I into the Hamiltonian and taking a functional derivative 



dyianlf)=Tj2 



(A.4) 



One can immediately calculate the derivative 



,5Sf(iw™) _ 1 



T^L - Unl 



(A.5) 



from the solution of the single-impurity problem. After substituting this result into Eq. (A. 3), we find that 



ma h 

m a ^"^ I fe fe' 



k k' ) 

- ^ E E - U^'f + - n))] 1 ^ Gr (io;™) (eL|-Q + 4+i-Q) 



The summations over momentum can be explicitly performed as follows: 
1 ^Gf (ia;„)efc_|_Q = X' [ZmF^{Z^) - l] ; 



1 ^Gr(ia;„)eL|_Q(efe_, - e^) = ^(1 - X) [Z^F^{Z^) - l] 



3 2 - - "2 - - ''*^ 



where 



Finally, the sum rule (first moment of the response function) contains two contributions 

I = Ik + In- 

The first contribution comes from the kinetic energy term 

Ik = 2(1 -X) J (iu;/(a;)Im|^ [Z{w)F^[Z{u>)] - l] 



-X 



/2 



*21 



-t*^ [Z{uj)F^[Z{uj)] - 1] - Z\lo) [Z{u:)F^[Z{uj)] - l] + 



and the second one comes from the potential-energy term 

+00 



f / *2 

Iu= 1 ckjf{io)ImJ2\ [S"(a;) - Un-}] ^—{1 - X'^)G%cj) + X'^Z-{u) [Z{u;)FM'^)] ' l] 

2 [^(u;)Fo,[Z(a;)]-l]'^ 



G'»(a;) 
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(A.6) 



(A.7) 



(A.8) 
(A.9) 

(A.10) 



(A.11) 



+ Uin} - n)) —(1 - X'2)G""(w) + X'^Z'^iu) [Z{uj)F^[Z{u;)] - l] 



Using the identities Z{uj)F^[Z{uj)\ — 1 = term contribution in the final form of Eq. (68), where 

\^{uj)G'^'^{lo) = X^{uj)G^^{uj) and [G'"'{lj)]''^ = we use the fact that in equiHbrium /x^ = /if and 
Z°-{u!) - X°-{uj), we can rewrite the potential-energy Z'^(a;) - Z'^(w) = I]^(a;) - E'^(w). 
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